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Abstract 

We study the four-dimensional SU(3) gauge model with a funda- 
mental and an adjoint plaquette term in the action. We investigate 
whether corrections to scaling can be reduced by using a negative value 
of the adjoint coupling. To this end, we have studied the finite tem- 
perature phase transition, the static potential and the mass of the O^"*" 
glueball. In order to compute these quantities we have implemented 
variance reduced estimators that have been proposed recently. Cor- 
rections to scaling are analysed in dimensionless combinations such as 
Tel \fo and mo++ jT^- We find that indeed the lattice artefacts in e.g. 
mo++ jTc can be reduced considerably compared with the pure Wilson 
(fundamental) gauge action at the same lattice spacing. 



1 Introduction 



Due to the enormous effort required by lattice QCD simulations, one is re- 
stricted to rather coarse lattice spacings. Therefore it is important to choose 
the lattice action such that already at a coarse lattice spacing the lattice 
artefacts are small. In the past, various proposals had been made to this 
end. 

One approach are so called perfect actions (for recent work on this subject 
see e.g. ref. Here one tries to avoid lattice artefacts completely, at the 
price of, in principle, an infinite number of terms in the action. In order to 
make this approach work practically, the set of terms has to be truncated to 
a finite, still large number of terms. The error introduced by this truncation 
is hard to estimate a priori. 

Alternatively, Symanzik ^ has proposed a scheme that allows to elim- 
inate lattice artefacts systematically, order by order in the lattice spacing 
and the coupling constant. For the SU(3) lattice gauge theory this has been 
worked out by Liischer and Weisz [2] up to the 1-loop level in perturbation 
theory to eliminate O(a^) corrections. 

Since we work at rather large values of the gauge-coupling it is not clear 
a priori, whether the 1-loop improvement of ref. is of any help. Only 
studies (see e.g. ref. of the scaling behaviour of various quantities can 
clarify this question. 

Here, we pursue a rather ad hoc approach. It is well established that the 
pure SU(3) gauge theory in four dimensions has a line of first order phase 
transitions with an end-point at 

(/?y,/5„) = (4.00(7), 2.06(8)) . (1) 

For the precise definition of /?/ and Pa see eq. ((Zj) below. At such transitions, 
lattice artefacts might completely disguise the continuum physics. Here, in 
particular, the mass mo++ of the lightest glueball 0"*"+ goes to zero as the end- 
point, eq. (H}, is approached 0. As a relic of this behaviour, the estimate 
of mo++/Tc (Tc is the location of the finite temperature phase transition) at 
/?a = 0, for 5.5 < (3f < 6.0 is much smaller than the continuum result. Here, 
we study whether the scaling behaviour can be improved by staying apart 
from the transition line by choosing negative values of f3a- Recently one of 
us Pj has studied gauge actions with 2x1- Wilson loops in addition to the 
plaquette. These actions are motivated by the RG-group and the Symanzik 
improvement programme. Indeed, for the so called Iwasaki action 7j an 
improved scaling behaviour was observed. Here we investigate whether by 
adding the adjoint part to the action similar benefits can be achieved, while 
keeping the action as local as possible; i.e. using plaquette terms only. For 
the mixed fundamental- adjoint action, it is straightforward to construct a 
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transfer matrix along the lines of ref. ^7j. However, for our choices of jSa-, it 
is not positive. For a detailed discussion see the appendix. 

An other important question (that we do not study here) is, whether 
dislocations, i.e. topological objects that are pure lattice artefacts, can be 
suppressed by a modification of the standard Wilson gauge action. 

The outline of our paper is as follows: First we give the definition of 
the model and the observables that we have studied. Next we discuss the 
Monte Carlo algorithm that was used to generate the gauge-configurations. 
In section 4 we study the finite temperature phase transition. In section 5 we 
extract the string tension a from the Polyakov loop correlation function. To 
this end we have implemented a variant of the variance reducing algorithm 
proposed by Liischer and Weisz jH]. Next, in section 6, we determine the 
mass mo++ of the lightest glueball O"*"^. Also here, we have employed a new 
algorithm [H], to reduce the variance of correlation functions. Based on these 
results, we compute the dimensionless quantities Td y/a and mo++/Tc. 



2 The model 

We consider a hypercubical lattice with the linear extension aNg in the spacial 
directions and aNt in the temporal direction, a is the lattice spacing. Ns 
and Nt are integer numbers. In all directions, periodic boundary conditions 
are applied. 

The action describing the pure SU(N) lattice gauge theory containing 
plaquette terms only, can be written in the general form 



s = E 



-ReTi aUp 



(2) 



where the sum over a indicates the sum over all representations of SU(N), 
Pa is the coupling associated with each representation, d{a) is the dimension 
of the representation, Tr^ is the trace in the given representation and Up is 
the path ordered product of the link variables along an elementary plaquette 
P: 

Up{x,n,u) = U^^^U^+^^Ml+u,f,Ul^ , (3) 

where x labels the sites of the lattice and /t is a unit vector in the /i-direction. 
In particular, we will consider the mixed fundamental- adjoint action 



S 



p 



-ReTTfUp 



+ Pa E 
P 



iV2 



-ReTr,[/p 



Using the identity 



Tr„t/ = Tr.f/^Tr.f/ - 1 



(4) 



(5) 
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and defining 
one obtains 

s = PfT. 



m - 1 



Pa , 



1 - —ReTTfUp 



/5a E 



1 



1 - — Tr/f/]>Tr/f/p 



(6) 



(7) 



which is the form that we will use in the following for SU(3). In the naive 
continuum limit the bare coupling qq for the action eq. (|7j) is given by 
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7o 



I3w = l3f + 2f3a 



(8) 



where we have introduced, following the literature, the notion of an "equiv- 
alent Wilson coupling" Pw 



2.1 Observables 

In the following we define the basic observables that we have studied. The 
fundamental plaquette is given by 



E 



f 



1 

QV 



E (Tr/f/p) , 



(9) 



X,fl,U>fl 



where V = N^Nf is the number of points of the lattice. 
The adjoint plaquette can be defined by 



(10) 



To study the finite temperature phase transition and the potential of static 
quarks we consider the Polyakov loop 



P,^ = Tr^nf^- 



i,0 



(11) 



where Ug^t,^l is the link variable starting in x, t in the direction /i. In partic- 
ular, we study the histogram associated with the quantity 



S X 



(12) 
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3 Simulation algorithm 



We have simulated the model with a local Metropolis procedure. The pro- 
posal for a new value of a link variable is generated by a Cabbibo-Marinari 
(CM) heat-bath update of a single SU(2)-subgroup for the action 



p 



1 - l^ReTifUp 



(13) 



This proposal is then accepted with the probability 

A = min [1, exp(-5(f/') + 5o(f/') + S{U) - SoiU))] , (14) 

where U is the original gauge field and U' is the proposal. < /?/ is tuned 
such that the optimal acceptance rate is obtained. In addition, we have 
performed overrelaxation (OV) updates [12] that keep the fundamental part 
of the action constant. Here, we accepted the proposal with the probability 

A = mm[l,exp{-Sa{U') + Sa{U))] , (15) 

where 

SaiU) = /?a E [l - ]^Tr;[/],Tr;[/p] . (16) 

In both cases we applied, for a given link, sub-group updates for the 1-2, 1-3 
and 2-3 components in a sequence. 

Using these elementary link-updates we are sweeping over the lattice. A 
complete update cycle is given by 

• One Cabbibo-Marinari-Metropolis sweep 

• M overrelaxation-Metropolis sweeps. 

We have implemented the algorithm in C as well as in Fortran. On the 
one hand, this is a good check for the correctness of the implementation and 
on the other, for the determination of the glueball mass, we intended to use 
part the Fortran code that was written for the study reported in ref. jlj. 
We have used the random number generator discussed in ref. For both 
implementations, we find that the time needed to update a link variable for 
the mixed action is roughly twice the time needed for the Wilson gauge action 
(i.e. (3a = 0). With the C-program, we need 1.4 x lO'^s (CM-Metropolis) 
and 0.9 x 10~^s (OV-Metropolis) for the update of a single link variable on 
a Pentium 4 PC running at 1.7 GHz. 
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Table 1: Acceptance rate as a function of jS'j. We have simulated a 4 x 16^ 
lattice at Pf = 9.25 and Pa = —4.0. An update-cycle consists of one CM- 
Metropolis sweep followed by 5 OV-Metropolis sweeps. In each, case, we 
performed 20000 cycles with 1000 discarded for thermalisation. 





Pace, heat 


p 

ace, over 


Ef 




1^1 


4.0 
5.0 
6.0 
7.0 
9.25 


0.76799(5) 
0.83647(4) 
0.82884(1) 
0.77845(2) 
0.64654(2) 


0.80855(4) 
0.80861(4) 
0.80854(3) 
0.80856(3) 
0.80855(3) 


1.6146(2) 
1.6148(2) 
1.6145(1) 
1.6145(1) 
1.6145(1) 


2.9472(5) 
2.9478(5) 
2.9469(4) 
2.9471(4) 
2.9470(4) 


0.056(5) 
0.062(5) 
0.055(4) 
0.056(4) 
0.054(4) 



3.1 Comparison with the Uterature 

As a test of the correctness of the program we tried to reproduce the values 
for the fundamental and adjoint plaquette given in fig. 2 and Table 1 of ref. 

From fig. 2 one reads off that Ef 1.35 for /3/ = 3.8, Pa = 2.25 and 
Ef ^ 1.88 for Pf = 4.0, pa = 2.25. Simulating an 8 x 12^ lattice we find 
Ef = 1.3465(12) and Ef = 1.8803(4), respectively. In both cases, 700 update 
cycles each with one CM-Metropolis sweep followed by M=5 OV-Metropolis 
sweeps were performed. We have discarded 200 and 100 cycles, respectively. 

In Table 1 of ref. the authors give values for the jump of the funda- 
mental and adjoint plaquette at the first order bulk phase transition. Among 
other values, they give AEf = 0.656(2) and AEa = 0.464(2) at pf = 3.27, 

Pa = 3.0. 

We computed these values on a 8 x 12^ lattice. We have started one 
simulation with an ordered and one with a disordered configuration. Within 
the simulation we did not observe tunneling between the phases. Hence we 
computed AEf and AEa as the difference of Ef and Ea obtained from the 
run with ordered and the run with disordered start. Our result is AEf = 
1.8863(7) - 1.2287(10) = 0.6576(12) and 3AEa = 3.947(2) - 2.550(2) = 
1.397(3) in perfect agreement with ref. [S]. Both simulations had 400 cycles 
with 50 cycles discarded for equilibration. 

3.2 Tuning the algorithm 

To this end, we have simulated a 4 x 16^ lattice at Pf = 9.25 and Pa = —4.0. 
In Table ^ we give the acceptance rate of the OV-Metropolis and the CM- 
Metropolis step as a function of Pf. We see that in the case of the CM- 
Metropolis the acceptance rate depends very weakly on P'f. The optimal 
value is reached for P'r around 5.0 to 6.0. The acceptance rate for the OV- 
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Metropolis is the same for all runs, as it should. With more than 80% it is 
still reasonable. 



4 The finite temperature phase transition 

The pure Yang-Mills theory undergoes a first order phase transition at a 
finite temperature Tc |14| IT5] . In the high temperature phase, the system 
is disordered and in the thermodynamic limit, the expectation value of the 
Polyakov loop vanishes. On the other hand, in the low temperature phase, 
the Z3 centre symmetry is broken. 

For a lattice with the extension of Nf lattice spacings in the time direction, 
in the limit Ns ^ 00, the deconfinement temperature is given by 



where {(3f,Pa}c indicates the critical coupling. In our numerical study, we 
could not compute the complete critical curve {Pf,l3a}c for a given value of 
Nt, but instead, we have only determined for the fixed values = 0, — 2 
and —4 of the adjoint coupling constant. 

In order to determine jSj^c, we have studied the probability distribution 
(for the definition of Q see eq. ()12p ). At the transition point, for 
sufficiently large lattice sizes, the histogram of the order parameter has a 
double peak structure. The weight of each of the phases should be the same 
at the transition point (see ref. [E])- One should note that, due to the Z3 
centre symmetry, the ordered phase is threefold degenerate. On the finite 
lattice we need some rule, how to assign given configurations to a particular 
phase. Here, we assign configurations with < Omin to the disordered 
phase, and configurations with \ > Omin to the ordered phase, where Omin 
is the minimum of p(|f2|) between the two peaks. As the lattice size increases, 
the separation of the peaks becomes sharper and therefore, the ambiguity of 
this assignment vanishes. 

In particular, we have computed the weight of the disordered phase as 
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(17) 




(18) 



and the weight of the three ordered phases as 




(19) 



The estimate of the critical (3f^c is given by 

Porderi^P f,c) g 

Pdis{Pf,c) 



(20) 
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where the factor 3 stems from the threefold degeneracy of the ordered phase. 
Following ref. ^Hl; the estimates for (3f^c should converge exponentially fast 
as Ns is increasing. 

In the Monte Carlo simulation, probabilities such as eq.s ()18pi9|l can be 
obtained in the following way: the probability that is in the interval 
[Oi,02] is estimated by 

1 ^max 

P[0,,0,]^ _ E X[o„o.](|f^|«) , (21) 

where X[0i,02] is the characteristic function of the interval [Oi,02]. is 
the value of for the i^^ configuration, N^ax is the total number of con- 
figurations that has been generated and Ndisc the number of configurations 
that have been discarded for equilibration. I.e., we just count how frequently 
the observable falls into the desired interval. 

Using reweighting we can easily obtain the histogram of |f2| for all values 
oi (3f m the neighborhood of where the simulation has been performed 
(i.e. the configurations are generated with the Boltzmann weight for I3f^s)- 
The standard reweighting method gives in our case: 

PfOi O^Md ) ~ ^£"^1.+! XtOi.02](l^l^^^) exp(-[/?/ - Y.x,,,u>, ReTr;[/ff ) 
' ^ ^ EiIX...+iexp(-[/^/-/^/JS-,/.,^>/.ReTr;f/«) 

(22) 

If (3f^s is close enough to /5/,c, i-e. a double peak can be seen at we use 
eq. (j22j) in combination with the intersection method to find the solution of 
eq. (EDD- 

Mostly, we have started our search for (3f^c on lattices with A^^ = 3A^t. 
Here it is helpful, to have a good first guess for [3f^c that can be used as a 
first simulation point [3f^s- To this end, see our discussion on hues of constant 
physics in subsection 14.31 

The result for (3f^c obtained with Ng = 3Nt was then used as a first guess 
Pf^s for the simulation of the next larger lattice size and so on. For lattice 
sizes Ng < QNt we have used the Monte Carlo algorithm as described in the 
previous section for the simulation. For A^^ = 6Nt, where our final results 
for jSf^c are taken from, we have employed the multicanonical method on 
top of it. A discussion is given in the following subsection. 



4.1 Enhancing the TunneUng Rate 

As the lattice size increases, the separation of the phases becomes more 
pronounced. While this allows for an unambiguous separation of the phases 
and hence of the transition temperature, it has adverse effects on the Monte 
Carlo simulation itself. Since in each elementary step of the update, the 
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configuration is only changed by a little, i.e. just at a single link, going form 
an ordered phase to the disordered, we have to pass through configurations 
that are mixtures of the phases. However, if such configurations are strongly 
suppressed, the Monte Carlo time to go from one phase to the other will 
be very large. In order to overcome this problem, it has been proposed not 
to generate the configurations with their Boltzmann-weight but rather with 
some modified one. In order to enhance the rate of tunneling events from 
one phase to the other we have simulated with a multicanonical ensemble 
[T7] . I.e. we have generated the configurations with a weight proportional to 

B[U]=exp{-S[U])W{\n\) , (23) 

where the modification factor only depends on the modulus of the 

sum over the Polyakov loops. 

We sweep over the whole lattice with the CM-Metropolis and the OV- 
Metropolis as described above. In order to fulfil detailed balance, with a 
probability of 50% the update sweep is performed in exactly reversed order. 
After these sweeps, an accept/reject step is performed with the acceptance 

A = mm[l,W{\Q'\)/W{\Q\)\ , (24) 

where is the value of for the proposed configuration. For the modifi- 
cation factor 14^(|f2|), we have used the ansatz 



W{ 


n\) 


= 1 




for 






\n\ 


< Sl 


W{ 


Q\) 


= i + m{\n\- 


Sl)/(S2 - Si) 


for 


Sl 


< 


\Q\ 


< S2 


W{ 


Q\) 


= 1 + m 




for 


S2 


< 


\Q\ 


< S3 


W{ 


Q\) 


= 1 + m (s4 — 


\n\)/is,-s3) 


for 


S3 


< 


\Q\ 


< S4 


W{ 


Q\) 


= 1 




for 


S4 


< 


\Q\ 





(25) 



The parameters si, S2, S3, S4, S5 and m of this ansatz were essentially deter- 
mined by trail and error. We have used this method only for our largest 
value of A^^: Ns = QNt for a given temporal extension Nt. 

4.2 Numerical results for /3f^c 

Let us discuss in detail the example of a 6x36^ lattice simulated at Pf = 7.803 
and f3a = —2.0. For this simulation, we have chosen the parameters of W as 
Sl = 0.03, S2 = 0.06, S3 = 0.08, S4 = 0.15, and m = 4.0. 

The simulation consists of 48000 update-cycles with one CM-Metropolis 
and five OV-Metropolis sweeps each. The acceptance rate for the modified 
weight eq. (f^^ was about 93%. In fig. ^the evolution of in Monte Carlo 
time is shown. We see that the system indeed tunnels quite frequently from 
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Table 2: Summary of our results for the finite temperature pliase transition 
obtained for tlie fundamental-adjoint gauge action. We give the value of (3f^c 
for fixed values of Nt and (3a- For a discussion of the numbers see the text. 



Nt; Pa 


0.0 


-2.0 


-4.0 


2 


5.0948(6) 


6.4475(6) 


7.8477(6) 


3 


5.5420(3) 


7.1603(3) 


8.8357(4) 


4 


5.6926(2) 


7.4433(3) 


9.2552(6) 


6 




7.8056(5) 


9.7748(11) 



Table 3: Results for the finite temperature phase transition obtained by 
other authors for the standard Wilson gauge action, i.e. = 0. 







ref. 


2 


5.0933(7) 


m 


3 






4 


5.6927(4) 


CHI 


4 


5.69254(24) 




4 


5.6925(2) 


EDI 


6 


5.89405(51) 




6 


5.8941(5) 


M 


8 


6.0624(9)(3) 




12 


6.3380(13)(10) 


20 ,[2T1 



the disordered phase into one of the ordered phases. For equilibration, we 
have discarded the first 2000 update-cycles. 

In fig. |21we give the histogram for exp{—S[U]) W{\Q\). Indeed, with our 
choice of parameters for iy(|f2|), no deep minimum, separating the peaks, 
is visible. Fig. El gives the Boltzmann-distribution for /?/ = 7.803 and 
Pa = —2.0 obtained from the simulation. Now we see a clear minimum that 
separates the peaks for the disordered phase and the three ordered phases. 
Still the weight of the two peaks is roughly the same. Therefore we per- 
formed a reweighting following eq. (j22I); such that eq. (j^UI) is satisfied. Under 
reweighting, the position of the minimum Omin slightly shifts. To get a con- 
sistent result, we therefore repeat the analysis with the new value of Omin- 
We found that the procedure quickly converges, such that we get an estimate 
of Pf^c with the proper Omin- The analysis of the statistical error of Pj^c is 
done with the standard Jacknife procedure. 

Our results are summarized in Table El For comparison, we have collected 
results for Pa = from the literature. In contrast to our study, these results 
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50000 



Monte Carlo Time 



Figure 1: Evolution of in Monte Carlo time. A unit in time is given by 
a cycle, consisting of one CM-Metropolis and 5 OV-Metropolis sweeps. We 
have simulated a 6 x 36^ lattice, at /?/ = 7.803, /3a = —2.0, with the modified 
weight that is specified in the text. 
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0.1 
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Figure 2: Histogram for the modified probability distribution for a 6 x 36^ 
lattice at /3/ = 7.803, /3a = —2.0. Details are given in the text. 
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Figure 3: Histogram for the Boltzmann distribution 6 x 36^ lattice at fif 
7.803, j3a = —2.0. Details are given in the text. 
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Figure 4: Histogram for the proper probability distribution 6 x 36'^ lattice, 
simulated at /?/ = 7.803, Pa = —2.0, reweighted to /3/ = 7.8056. Details are 
given in the text. 

where obtained from the peak of the susceptibility. In the two cases, where 
we have mapping parameters, our results are consistent with those of the 
literature. For Pa = 0, Nt = 6 we performed no own simulation, but used in 
the following the resuh Pf = 5.89405(51) of ref. 

4.3 Lines of constant physics 

At first order in perturbation theory the tree-level relation (jH)) is modified as 

Up to higher order perturbative corrections and non-perturbative contribu- 
tions, lines of constant Pw should represent lines of constant physics. As a 
check of this relation, we compare in fig. Elthe values of the critical couplings 
for Nt = 6, i.e. a ~ 0.11 fm, * with the one-loop prediction eq. (jSH)), where we 
have set Pw = /5/,c|a=o.iifm, /3a=o- For Pa > we adopt the results obtained 
in 0. 

We see that the perturbative prediction fails in describing the line of constant 
physics. In ref. jSj it was observed that the prediction for positive Pa can 
be considerably improved by using a tadpole improved perturbative formula. 

*In order to express the lattice spacing in physical units, we used the relation TcTq = 
0.7498(50) 4 and ro ~ 0.5 fm 
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Figure 5: Line of constant physics (solid line) as predicted by eq. and 
the deconfinement transition line for a = 0.11 fm (open circles). 
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Figure 6: Deconfinement transition lines for a = 0.11, 0.17, 0.33 fm and the 
interpolation formulas eq. (solid lines). For /^^ > the data are taken 
from ref. [S]. 

For negative however, the perturbative prediction seems to be completely 
unreliable. Hence we made no attempt to study tadpole improvement here. 

However, it turns out that the deconfinement transition lines can be 
parametrised quite well by a quadratic interpolation. Fitting the full range 
of available Pa (including the results of ref. jH] for Pa > 0) yields: 

a = 0.33fm:/3f = 5.0948 - 0.6645/3^ + 0.0059/5^ (27) 
a = 0.17fm:/5f = 5.6930 - 0.8586/5^ + 0.0080/5^ 
a = 0.11 fm : /3f = 5.8951 - 0.9391/3^ + 0.0078/3^ . 

Fig. ini shows the deconfinement transition lines for a = 0.11, 0.17, 0.33 fm for 
a range of positive and negative Pa- For a = 0.33 fm, the interpolation for- 
mula reproduces the numerical values for Pa < within the statistical errors. 
For a = 0.17, 0.11 fm and /3a < the quadratic interpolation reproduced the 
numerical data only within two standard deviations. 
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5 The static potential 



We have extracted the static potential from Polyakov loop correlation func- 
tions: 

a\/(r) = -i-[ln(P(xrP(y)) + 6] , (28) 

where y = x + rl. e is the correction due to excited states in the string 
spectrum. It vanishes exponentially as A^^ — * oo. In the free bosonic string 
approximation one gets [2S1 123 1211 1211 

oo 

e = 2 ^ ln(l - expf-vmaiVt/r]) . (29) 

n=l 

Note that here, in contrast to the previous section, we are considering sys- 
tems with a temporal extension aNt » 1 /Tc to eliminate finite temperature 
corrections e as much as possible. In the simulations reported below, we have 
chosen Nt = 6/{aTc). It is straightforward to check that for this choice, at 
our level of numerical precision, the corrections eq. pUj) can be safely ignored. 

We have computed the Polyakov loop correlation function with a variant 
of the algorithm that was recently proposed by Liischer and Weisz jS]. Details 
of the algorithm are given below. 

In this study, we consider rather coarse lattice spacings. Therefore the 
Sommer scale |221 is intrinsically affected by large systematic errors. On 
the other hand, we have computed the static potential up to rather large 
physical distances r. Therefore we decided to compute the string tension a^a 
rather than ro/a. We evaluated the string tension using the ansatz 

nr)=ar + /.-^(^l + ^) (30) 

for the static quark potential, with b = 0.04fm where contributions of 
the order 0{r~^) are neglected. Details of our numerical analysis are given 
in the subsection 15.21 



5.1 Variant of the Liischer and Weisz method to com- 
pute the Polyakov loop correlation function 

In order to reach large values of Nt we have implemented a variant of the 
recent proposal of Liischer and Weisz jH] . In addition to the factorisation in 
temporal direction, we also use a factorisation in the spacial directions ^ 

'I'ln ref. |27] the Polyakov loop correlation function was measured with a spatial de- 
composition only. The model studied in this ref. contains scalar fields in addition to the 
gauge field 
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In temporal direction, we have used only one level of the factorisation. 
To this end, the lattice is divided in temporal direction into Nt/Ni layers of 
the thickness Ni. Note that Liischer and Weisz [Hj also have used only one 
level of the factorisation in most of their numerical studies. 

In addition, we have used a factorisation in the spacial directions. To 
this end, we have divided the lattice in blocks of the size x Ni. Within 
these blocks we consider only a subset of all possible Polyakov loops. See 
the two-dimensional sketch fig. [7|for an illustration. The idea of this choice 
is to use, for a given distance of the loops, only those loops that have a 
maximal distance from the boundaries of the blocks. In particular, for an 
even distance r/a between the loops, we took the two loops with distance 
(r/a)/2 from the common boundary. For an odd distance r/a between the 
loops, we took the distances (r/a+ l)/2 and {r/a — l)/2 from the boundary 
between the blocks. 

We have now a two-fold hierarchy of the algorithm. At the lowest level, 
we update at fixed boundaries of the blocks and fixed boundaries between 
the temporal slices. This step provides us with variance reduced segments of 
the Polyakov loops 

1 Muock tQ + Ni-1 

p{s,to) = --- E n f4:i,o , (31) 



Mblock 



where Muock is the number of updates that have been performed with fixed 
boundaries of the block and i labels the configurations that have been gener- 
ated this way. These variance reduced segments of the Polyakov loop could be 
viewed as a generalization of the multi-hit method for the variance reduction 
of a single link variable that has been applied in ref. jH]. 

Two of these variance reduced segments with the same to from neighbour- 
ing spatial blocks are now used to construct the complex 9x9 matrices of 
eq. (3.2) of ref . ^: 

T(f , y, to)a,/3,7,5 = ^(^> ^o)a,/3^(y, ^^0)7,5 • (32) 

For fixed boundaries between the temporal layers, we perform a certain 
number of update sweeps before we repeat again the calculation of the vari- 
ance reduced segments eq. (jSU- This procedure is performed Miayer times 
and the matrix T is averaged over these Miayer instances: 

J{x,y,to)a,0,'r,5 = jii T(f,y,to)S,T,,5 • (33) 

^^J-layer j—\ 

As in ref. |j8j, the Polyakov loop correlation function is now computed as 

1 Mmeas [ t=Nt/Ni-l _ ~\ 

(PixyPiy))^-—- Y: n nx,y,tNiY''n , (34) 

^^J-meas f—n 
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where Mmeas is the number of complete measurement cycles. It is understood 
that measurements start after equilibration of the system. 

In order to further improve the measurement of the segments of the 
Polyakov loop eq. (jHlll we have performed a multi-hit update of the single link 
variables. In addition, we have updated the links close to the centre of the 
block more frequently than those close to the boundary. To this end, we have 
set up a sequence of blocks inside the block having the sizes (6 — 2)^, (6 — 4)^, 
... . For each cycle (that consists of 2 OV-Metropolis and 2 CM-Metropolis 
sweeps over a block of the spatial size {b — 2m)^) we have performed n such 
cycles for the next smaller block (6 — 2(m — 1))^. I.e. for each cycle of the full 
block (b^), additional n, n"^, n^... cycles are performed for the blocks of the 
spatial size (6 — 2)^, (6 — 4)^, (6 — 6)^, ... inside the block of the spatial size 
b^. In most of our simulations we have chosen n = 3. This way, of course, 
we can not gain exponentially, but we can fight, to some extend, the factor 
that we loose by the reduced number of Polyakov loops that we consider. 

Compared with the original approach of Liischer and Weisz, we gain the 
factorisation in the spacial directions. However, we lose a lot of copies of 
the correlator. Instead of 3 x only 3 x (Ng/b)^ remain. To compensate 
partially for this fact, we update more frequently the links in the centre of 
the block than those at the boundary. A positive side-effect of the reduced 
number of copies of the Polyakov loop is that the memory requirements are 
drastically reduced compared with the original proposal. 

Most of the CPU-time is spent with the block-boundaries fixed. This 
would allow for a rather trivial parallelisation on {Nt/Ni){Ns/bY nodes. This 
makes the algorithm ideally suited for a PC-cluster equipped with a moder- 
ately fast interconnect like Gigabit-ethernet. 

The obvious disadvantage of our variant of the algorithm is that even 
more parameters have to be tuned as in the original one of Liischer and 
Weisz. In our case, it is almost unavoidable to make just ad hoc choices for 
some of the parameters. Finally, for the distances r/a that we have reached, 
our method still seems to be slightly less efficient than the one of ref. jH]. 
Here is difficult to give a definite comparison, since we did not simulate at 
exactly the same parameters as ref. [S]. 

5.2 Numerical results for the string tension 

To eliminate the constant fi in eq. (|3(J|) we consider the so called force 




(35) 



On the lattice, we compute the force from the potential either as 



aF{r - a/2) = V{r) - V{r - a) 



(36) 
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Figure 7: Two dimensional skctcli of tiie block. Two of the spatial direction 
are show. The remaining spatial direction and the the temporal direction 
are perpendicular to the paper. Only the Polyakov loops that are indicated 
with black discs are computed. 
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Table 4: Parameters of our simulations to compute the static quark po- 
tential. In the first two columns we give the coupling constants. The third 
column contains the lattice sizes. In the fourth column the range of the dis- 
tance r is shown where the static potential has been computed. Finally we 
give the number of measurements Mmeas, as defined in the previous section. 



Pa 


Pf 


Ns, 


Nt 


r/a 







5.0848 


16, 


12 


2-6 


30 


-2 


6.4475 


16, 


12 


2-6 


26 


-4 


7.8477 


16, 


12 


2-6 


24 





5.5420 


20, 


18 


2-8 


17 


-2 


7.1603 


20, 


18 


2-8 


17 


-4 


8.8357 


20, 


18 


2-8 


21 





5.6926 


24, 


24 


2-10 


56 


-2 


7.4433 


18, 


24 


2-4 


76 


-4 


9.2564 


18, 


24 


2-4 


101 


-4 


9.254 


18, 


24 


2-4 


76 



or 

aF{ri) = V{r) - V{r - a) , (37) 

where r/ is the tree- level improved distance defined in ref. [22] . 

In our numerical analysis, we made no attempt to compute b, but instead 
have used the result b = 0.04fm of ref. jH]. In order to obtain the dimension- 
less quantity a~^b we have used TcVq = 0.7498(50) obtained in ref. |1] with 
To = 0.5fm. It follows 

a{{Pf,Pa}cy'b^0.06Nt , (38) 

where here a({/?/, Pa}c)Nt = l/Tc. As a result, eq. has no free parameters 
in addition to a. Hence, for each value of r/a we obtain an estimate for a^cr. 

We computed the static quark potential at the critical couplings that we 
have evaluated in the previous section for ^ = 2,3,4. The parameters of 
these simulations are given in Table |3 As thickness of the temporal layers 
we have used A^^^ = l/(aTc) throughout. For the other parameters, let us just 
detail a typical example: For = 0, Pf = 5.6926 we have used Miayer = 30 
and Mhiock = 40. 

In Table El we give details of our analysis of the force for Pa = 0, Pf = 
5.6926, where we have collected our most accurate data. First we have esti- 
mated the string tension as a{V{r) — V{r — a)), ignoring the Liischer term. 
We see that even for our largest values of r/a the estimate does not stabilize. 
This behaviour is clearly improved, when the Liischer term is used in the 
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Table 5: Analysis of the force at (3a = 0, (3f = 5.6926. In this Table we give 
results obtained with various ansatze. In the second column we give the most 
naive estimate for the string tension o?a: a{y{r) — V{r — a)). In the third, 
we use the ansatz J^Sj) in connection with eq. and 6 = 0. In the fourth 
column, the naive definition of the distance is replaced by the improved one 
eq. n^) . still 6 = 0. Finally, in the column 5 we use b = 0.04fm as discussed 
in the text. In the last column we report the statistical error, which is the 
same in all cases. 



r/a 


0- ^naive 






2 


stat. error 


3 


0.2211 


0.1793 


0.1707 


0.1600 


0.0003 


4 


0.1902 


0.1688 


0.1664 


0.1629 


0.0004 


5 


0.1783 


0.1654 


0.1645 


0.1630 


0.0004 


6 


0.1727 


0.1640 


0.1637 


0.1629 


0.0006 


7 


0.1692 


0.1630 


0.1628 


0.1623 


0.0006 


8 


0.1667 


0.1621 


0.1620 


0.1617 


0.0011 


9 


0.1655 


0.1618 


0.1618 


0.1616 


0.0025 



ansatz. Fortunately, the difference between the two definitions of the argu- 
ment of the force eq.s ()36I37|) is only minor. Finally, we included b = 0.04fm 
in the ansatz. As a result, the estimate of the string tension further stabilizes 
as the distance r/a is varied. In fact, starting from r/a = 4, all results are 
consistent within the statistical error. 

As our final result we quote a^a = 0.162(1), which is obtained from the 
ansatz with r/ and b = 0.04fm at the distance r/a = 7. Since the difference 
of our results with 6 = and b = 0.04fm at r/a = 7 is smaller than the 
statistical error, we are confident that the quoted error also covers possible 
systematic errors. 

We have extracted the final result for a^cr for the other values of the 
couplings in a similar fashion. In particular, we have always taken the result 
obtained with the tree-level improved distance rj and with b = 0.04fm. For 
l/(aTc) = 2 and 3, the final results for a^cr are taken from r/a = 4 and 
r/a = 6, respectively. In the case of l/(aTc) = 4 we give the final result for 
Pa = and Pf = 5.6926 that we obtained above. For Pa = —2 and Pa = —4 
we do not have data for such large distances. Therefore we took r/a = 4. 
The error that is quoted here includes a systematic error that we estimate 
based on our data for = and Pj = 5.6926. These results are summarized 
in Table ini 

M. Liischer j2Hl has provided us with numerical data for the force at 
Pf = 5.7, Pa = 0.0 obtained on a 24 x 18^ lattice |H]. Using the procedure 
discussed above, we extract a^a = 0.156(1) from these data, where the error 
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Table 6: Summary of our final results for the string tension a? a and the 
dimensionless ratio T^/ ^/a. For a detailed discussion see the text. 



l/(aTe) 


Pa 


/5/ 






2 





5.0948 


0.759(2) 


0.574(1) 




-2 


6.4475 


0.742(2) 


0.580(1) 




-4 


7.8477 


0.736(2) 


0.583(1) 


3 





5.5420 


0.319(2) 


0.590(2) 




-2 


7.1603 


0.308(2) 


0.601(2) 




-4 


8.8357 


0.306(2) 


0.603(2) 


4 





5.6926 


0.162(1) 


0.621(2) 




-2 


7.4433 


0.160(1) 


0.625(2) 




-4 


9.2540 


0.160(1) 






-4 


9.2564 


0.159(1) 






-4 


9.2552 




0.626(2) 



mainly covers possible systematic errors. (The statistical error of the force 
at r/a = 7 is only 1.5 x 10""^). Extrapolating our data for Pf = 5.5420, 

= and Pf = 5.6926, Pa = 0, using the ansatz In(a^cr) = c + dp, we get 
a^cr = 0.1567(10) for Pf = 5.7, Pa = 0, which is consistent with the result 
that we have extracted from the numerical data of ref. |H1 128j. 

In Table ini we also give the results for Tc/i/a, which are plotted in fig. IHl 
together with other values obtained for the Wilson action in [21]. The error 
of Tc/ \fo is dominated by the error of a^o. 

We see that for Pa = —2 and —4 the estimate for Tc/ \fa is closer to the 
continuum limit than for Pa = 0. However, the difference between l/(aTc) = 
3 and l/(aTc) = 4 is larger than that for the different values of Pa at fixed 
l/{aTc). Since already for l/{aTc) = 4 there is only a minor difference in the 
estimates for \fo from /5a = and Pa = —2,-4, we decided to skip the 
determination of a^a for l/{aTc) = 6. 

For comparison, we have summarised in Table [7| continuum limit results 
for Tc/ \fa given in the literature that are obtained with various actions. In 
particular, the result of ref.s |^ and [HI] are not compatible within the given 
error bars. Given that our smallest lattice spacing is a ~ 0.17 fm, we did 
not extrapolate our data to the continuum limit. However, looking at fig. |H1 
our data seem to be more compatible with a continuum result Tc/ \fo ~ 0.66 
than 0.63. 
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Figure 8: We have plotted our results for Tc/ y/a as a function of where 
Nt = l/{aTc). In addition we give the results obtained in ref. |2I] for Pa = 
at smaller lattice spacings. 



Table 7: Results for the continuum limit of T^/ \fo obtained with various 
lattice actions. 



action 




Wilson [21] 
Symanzik imp. [21J 
DBW2 |2i 
Iwasaki [201 

1-loop tadpole impr. |3T] 


0.630(5) 

0.634(8) 

0.627(12) 

0.651(12) 

0.659(8) 
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Table 8: Parameters of the glueball computation. For definitions see tlie 
text. 





nf 




Ns 


Nt 


ni 


Nor 


Nsub 


f-sub 


hi 


Nmeas 


5.0948 







4 


6 


1,2,3,4 


2 


320 


4 


6 


9931 


0.44/^0 


4.0 


— Z.U 


4 


p. 




1 O Q /I 

1,2,3,4 


Z 




Z 





1 nonn 
iyoUU 


7.8477 


4.5 


-4.0 


4 


6 


1,2,3,4 


2 


160 


2 


6 


17800 


5.5420 







6 


8 


2,4,6,8 


4 


160 


4 


8 


4894 


7.1603 


5.0 


-2.0 


6 


8 


2,4,6,8 


4 


80 


2 


8 


8500 


8.8357 


5.0 


-4.0 


6 


8 


2,4,6,8 


4 


80 


2 


8 


7140 


5.6926 







8 


12 


2,4,6,8 


5 


160 


4 


8 


4520 


7.4433 


5.0 


-2.0 


8 


12 


2,4,6,8 


5 


80 


2 


8 


2925 


9.2564 


5.0 


-4.0 


8 


12 


2,4,6,8 


5 


80 


2 


8 


3852 


5.89405 







12 


18 


3,6,9,12 


7 


300 


6 


10 


584 


7.8056 


5.0 


-2.0 


12 


18 


3,6,9,12 


7 


100 


2 


10 


672 


9.7748 


5.0 


-4.0 


12 


18 


3,6,9,12 


7 


100 


2 


10 


833 



6 Glueball masses 

We continued our investigation of the lattice artefacts of the mixed action 
by measuring the mass (at zero temperature) of the lightest glueball 0"'""'' 
at the critical values of [3f for [3a = 0,-2,-4 and l/(aTc) = 2,3,4,6 that 
we have determined in section 0J This allows us to study the scaling of 
the dimensionless quantity mQ++/Tc. The 0++ glueball mass is particularly 
interesting since it shows large lattice artefacts in the case of the Wilson 
action [j3a = 0). The fact that the 0"'"+ mass becomes very small at certain 
lattice spacings can be interpreted as the influence of the endpoint of the first 
order phase transition (see eq. (jT])), where mo++ vanishes jB]. If this picture 
is correct, we expect that by choosing a negative jSa, we move away from the 
endpoint and hence the lattice artefacts on mo++ should be reduced. 

In ref. [23 an action with 4 terms had been used: In addition to the 
plaquette in the fundamental representation, the plaquette in the adjoint 
representation, the represent ant ion of dimension 6 and a ltimes2 loop in the 
fundamental representation were considered. All couplings except that of the 
fundamental plaquette were negative. The authors argue that this way the 
end-point of line of first order phase transtions in the fundamental-adjoint 
can be avoided. Their final estimate for the mass of the lightest glueball is 
mo++ / sqrt{a) = 3.5(3), which corresponds to mo++ro = 4.07(35). 

On anisotropic lattices {at << a^) the authors of ref. ^| found that by 
introducing in the action a term ReTr?7p(a:, t, fi, z/) ReTiU p{vecx, t + at, fi,^) 
for space-like plaquettes with a negative coupling constant, the scaling be- 
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Figure 9: Wilson loops used for the evaluation of mo++. 

liaviour of the 0^"^ glueball mass can be improved. On top of this modification 
they jni] employed Symanzik-improvement |21 E] with a 2 x 1 Wilson loop 
term in the action. 

In order to extract the mass mo++ of the 0"*"+ glueball, we have computed 
the connected correlation function between spatial Wilson loops. To this 
end, we have chosen a basis of = 7 operators in the Af^ representation of 
the cubic group (among the 22 spatial Wilson loops up to lenght 8) plotted 
in fig. El which had shown the best signal-to-noise ratio in a previous study 

111 

In order to get a better overlap with light states, we applied the APE 
smearing procedure to the spatial links jHSj- For each operator, we used 
M = 4 different smearing levels. For each shape with d different orientations, 
we measured the observable 

M-3/2 d 

Oiit) = ^j:Y.T^rWrix,t) , (39) 

V« n=l 

where Wp is the spatial Wilson loop with smearing level /, in the orientation 
n. We then constructed the correlation matrix 

a,{t) = mt)om) - mt)){0M) , m 

where now the indices i,j run from 1 to x M = 28. 
6.1 Factorization and error reduction 

The exponential error reduction proposed by Liischer and Weisz jH| and 
adopted in our work to compute Polyakov loop correlation functions, can 
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be applied to a wider class of n-point functions. ^ In particular, this idea has 
been tested on the Wilson loop correlation functions to compute the scalar 
and tensor glueball masses jH]. 

Here we have implemented only one level of factorisation. To this end, 
we divided the lattice in temporal direction into two sub-lattices with an 
extension = Nt/2 each. While keeping the spatial links at the time slices 
to = 1 and ti = Nt/2 + 1 fixed ^, we performed Msub "sub-measurements" 

_ 1 J^su b 

O^it) = E . (41) 

JMsub n=l 

For each of these measurements, we perform Igub sub-updating cycles. Every 
updating consists of one Cabibbo-Marinari heatbath sweep and N^r overre- 
laxation sweeps. "Sub-updating" means that all link-variables except those 
at to and ti are updated. The total number of sub-update cycles with one 
set of fixed link variables at to aiid ti is then Ns^i, = Isub^sub- 

After these Ng^b sub- update cycles with fixed links at to and ti, we per- 
formed Igi global sweeps, with the same ratio Nor between overrelaxation 
and heatbath sweeps. For t even, we evaluated the 2-point "unsubtracted" 
correlation function as 

1 ^meas "1 

W— ^9 ^ A (f + t/2) O, (f - t/2) (42) 

JVmeas j=i ^ t'=l,Nt/2+l 

while for t > 2, t odd, we used 

1 ^meas 1 

Cr.it) ^ —— E I E 0,{t' + (t + l)/2) 0,(t' - (t - l)/2) 

^Vmeas j=i t'=l,Nt/2+l 

+ 0,(t' + (t - l)/2) Oj{t' - (t + l)/2) , (43) 

where Nmeas is the number of times eq. ()41|) is evaluated. An important 
point on the subtraction of the vacuum expectation value is that only the 
measurements included in Cfj{t) are taken into account. E.g. for t even, we 
subtract 

1 I Nmeas \ (Nmeas \ 

T7^ E A(t' + t/2) 0,(t'-t/2) . (44) 

In this way, the two terms are highly statistically correlated, and hence sta- 
tistical fluctuations cancel, when the difference is taken. For a more detailed 
discussion of this point see ref. 0. Alternatively one could consider a tem- 
poral derivative of the correlation function, as proposed in ref. |45j . 

■''A similar idea had been exploited in ref. |36| to compute fermion propagators. 
§In this section, we have set a = 1 to simplify the notation 
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The free parameter of this 2-level scheme is the number of sub-measurements 
Msub- Here, we did not try to tune Msub but rather followed the rule of ref. 

Kub^e^' , (45) 

where m is the mass of the state that should be measured and t is the time, 
where we intend to extract the mass from the exponential decay of the 2-point 
function. 

In fig. ^1 we analyse the effective mass evaluated with the usual algo- 
rithm and with the factorisation formula for /3/ = 5.6926, (3a = 0, up to 
t = 4. For this purpose we collected 16000 measurements performed with 
the usual method and 400 measurements obtained with one level of factori- 
sation, where each measurement is obtained with 40 sub-measurements. The 
computational effort for these two simulations is then roughly the same. The 
masses were extracted with the variational method discussed below. One can 
notice that for our choice of the algorithm, for t = 2 the standard method is 
still more efficient that the variance reducing one. In our final analysis, we 
did not make use of the factorisation method for time separations t = 0,1 
and we computed the correlation function in the usual way. For t > 3 we 
observe however a substantial error reduction, and actually this is the region 
where we extract the glueball mass. 

6.2 Analysis and results 

As a first step, we symmetrised the correlation matrix by replacing Cij{t) by 
1/2 [Cij{t) + Cji{t)]. The masses were extracted from the correlation func- 
tions by applying the usual variational method jSZlEHl- We solved the gen- 
eralised eigenvalue problem 

C{t) v^{t, to) = A„(t, to) C{to) v^{t, to) , (46) 

with a = 0, (A^ x M) — 1 and Aq > Ai > ... > A(ivxAf)-i- The indices i,j of 
the correlation matrix C{t) are now omitted. In our analysis, we have chosen 
to = throughout. We projected the correlation matrix on the ground state 
eigenvector computed for t = to + 1 

Wit) = v'^ito + 1, to) C(t) voito + 1, to) , (47) 

and evaluated the effective mass using 

rueffit) = log (W{t - 1) + ^W{t - 1)2 - W{Nt/2)^^ - (48) 
log (W{t) + sjw{t)^ - W{Nt/2)^^ , 
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Figure 10: The plot shows the effective mass of the 0"'""'" glueball at /?/ = 
5.6926, j3a = 0. In one case (squares) the correlation matrix has been cal- 
culated in the standard way. while in the other case (triangles) the variance 
reduction method that is discussed in the text has been apphed. Note that 
ior t — 1 also the triangle result is computed in the standard way. 
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Table 9: Results for the 0"*"+ glueball mass in lattice units. As our final result, 
we liave taken the effective mass evaluated at the distance t. 





Pa 


t 


amo++ 


mo++/Tc 


5.0948 





3 


2.237(88) 


4.47(18) 


bAAib 


-2.0 


2 


2.495(20) 


A r\r\r\ / a r\\ 

4.990(40) 


7.8477 


-4.0 


2 


2.645(25) 


5.290(50) 


5.5420 





3 


1.158(18) 


3.474(54) 


7.1603 


-2.0 


2 


1.414(14) 


4.242(42) 


8.8357 


-4.0 


2 


1.550(17) 


4.650(51) 


5.6926 





4 


0.967(16) 


3.868(64) 


7.4433 


-2.0 


3 


1.108(18) 


4.432(72) 


9.2564 


-4.0 


3 


1.193(18) 


4.772(72) 


5.89405 


4 





0.787(18) 


4.72(11) 


7.8056 


-2.0 


3 


0.839(18) 


5.03(11) 


9.7748 


-4.0 


3 


0.836(17) 


5.02(10) 



which takes into account the periodic boundary conditions in temporal di- 
rection. In fig. ^2we show as an example the effective mass for f3f = 5.6926, 
/3a = as a function of the distance t. We decided to extract the mass at 
t = 4, where the contributions of excited states should be smaller than the 
statistical errors. (Note that the mass of the first exited state in the Ai~^ 
channel is almost twice as heavy as the O"*"^ state.) Our final results are 
summarised in Table together with the value of t, where we extracted the 
masses. 

Our final results for mQ++/Tc are also reported in Table IHl Here the 
dominant error is the uncertainty on mo++. Note that we computed the 
glueball mass for j3f = 9.2564, /Jq = —4.0, which was a preliminary estimate 
for the critical coupling jSf^c and not the final result reported in Table |21 
However, we estimate that the shift in the O"*""*" glueball mass corresponding 
to the shift in /?/ is negligible with respect to the statistical error. Fig. IT^ 
shows mQ++/Tc as function of (aTc)^. 

For comparison, we give an estimate of the continuum limit based on the 
average of the following results given in the literature: 



mo++ro 


ref. 


4.35(11) 




4.33(10) 


P 


4.21(11)(4) 




4.23(22) 


[43^ 


4.30(6) 


average 



The computation presented in 02] has been performed with anisotropic lat- 
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Figure 11: Effective mass for the 0"*"+ state, at /?/ = 5.6926, = as a 
function of the distance t. 

tices. The resuhs of jlO] and [H] at a finite lattice spacing have been ex- 
pressed in units of tq in ref. |1H]- In ref. further results are reported 
for which the continuum extrapolation has not been performed. Results for 
mo++ro at finite lattice spacing obtained with the FP (fixed point) action are 
also present The error that we give for the average should not taken 

too serious, since it is not clear, to which extend the error of the individ- 
ual results is of systematic or statistic nature. In ref. [33] the authors plot 
their results for mo++ro obtained from anisotropic lattices with an improved 
action, as discussed above. They given no final result for the continuum 
limit. However, from the plot one reads off mo++ro ~ 4.0 with a quite small 
error; incompatible with the average of the literature given above by several 
standard deviations. 

Using the continuum limit relation [1] 

T,ro = 0.7498(50) (49) 

the average of the results from the literature can be converted to 

mo++/T,U=o = 5.73(9) . (50) 

We made no attempt to extract a continuum result from our data, since it is 
quite clear from fig. ^Jthat corrections beyond are large. At a ~ 0.1 Ifm 
we do observe a moderate reduction of the lattice artefacts by using < 
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Figure 12: mo++/Tc for /?a = (filled squares) Pa = —2 (open squares) and 
Pa = —4 (triangles) as function of (aTc)^. The circle gives the continuum 
results extracted from the literature. 



with respect to the usual Wilson action (/?„ = 0). For Pa = 0, the deviation 
from the continuum result of eq. (j^UI) amounts to ~ 18%, while for Pa = 
—2, —4 it slightly decreases to ~ 12%. 

At a ~ 0.17fm one observes discretization errors of ~ 40% for the Wilson 
action, while for the mixed action they amount to ~ 25% for Pa = —2 and 
~ 20% for Pa = -4. 

7 Summary and conclusions 

We investigated the SU(3) lattice gauge model with a pure gauge action 
that contains plaquette terms in the fundamental and adjoint representa- 
tion. In particular, we studied negative values of the adjoint coupling Pa- 
This choice is motivated by the presence of a first order phase transition line 
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in the {Pf,Pa) plane: one expects that moving towards negative Pa the pres- 
ence of the endpoint at {Pf,Pa) = (4.00(7), 2.06(8)) becomes less important 
and hence that scaling and/or topological properties could be improved with 
respect to the Wilson action. These features would be highly desirable in 
view of upcoming simulations for fermionic actions with exact lattice chiral 
symmetry and in general for unquenched computations. 

In section 4, we computed the critical coupling c of the finite tempera- 
ture deconfinement transition at l/{aTc) = 2, 3, 4, 6 and = 0, —2, —4 fixed. 
Since this measurement turned to be our most accurate, we have used Tc to 
set the scale. We find that hues of constant physics, as predicted by one- 
loop perturbation theory, completely fail to describe our numerical results 
for negative 

We then calculated the static quark potential from the Polyakov loop cor- 
relation funtion. To this end, we have implemented a variant of the algorithm 
recently proposed by Liischer and Weisz. In addition to the factorisation in 
temporal direction, we have employed a factorisation in the spacial direc- 
tions. This algorithm allowed us to compute the Polyakov loop correlation 
funtion up to r ~ 1.5fm. Due to these large distances we were able to extract 
the string tension a with little systematic errors (section 5). 

Studying the scaling behaviour of the quantity Tc/ y/a for the different 
Pa at our disposal, we did not observe a significant improvement at negative 
adjoint couplings in comparison to the Wilson case /3a = 0. The values 
obtained with negative f3a are a little closer to the continuum limit. 

In section 6 we computed the O"*"^ glueball mass for several lattice spacings 
and for = 0, —2, —4. Also for this computation we made use of the 
factorisation method to reduce the variance of the correlation function [Hj. 
Although the efficiency here is not as spectacular as for the Polyakov loop 
correlation function, we were able to obtain a good statistical accuracy up 
to distances roughly twice the ones reached with the standard method. 

It turns out that the mass mo++ of the lightest glueball is more sensitive 
to the variation of Pa- This had to be expected, since at the endpoint of 
the line of first order phase transitions the mass in lattice units is zero 
Therefore, in particular for mo++ large lattice artefacts should show up in 
the neighbourhood of the endpoint. We investigated the scaling behaviour of 
the dimensionless quantity mo++/Tc. Here indeed, we observed a significant 
reduction of the lattice artefacts for negative Pa- At a ~ 0.17fm, the lattice 
artefacts for Pa = are 40%, while for Pa = —4 they decrease to 20%. 

In view of future dynamical QCD simulations, it would be interesting to 
study the effect on dislocations and to investigate spectrum of the Wilson- 
Dirac matrix obtained with the mixed action at negative Pa- Finally one 
should note that the hybrid-Monte- Carlo (HMC) algorithm can be easily 
implemented for the mixed fundamental-adjoint action. 
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A Appendix: Is the transfer matrix positive? 

In ref. jUj the transfer matrix for lattice QCD with the Wilson (funda- 
mental) gauge action is constructed. It is straightforward to generalise this 
construction to the mixed fundamental /adjoint plaquette action. 

In ref. jlTj it is shown that the transfer matrix for the Wilson action is 
strictly positive if and only if 

JduJ dU'nU) exp [^^[Tr(f/-if/') + TT{U''U'y]j f{U') > (51) 

for all square integrable, nonvanishing functions / on the gauge group SU(3). 
This can be generalized as 

JdU J dU'f*{U) exp (^^[TrV^ + TiV^ + j^^^V Tr0^ f{U') > 

(52) 

where V = U'^U' . 

Following ref. |17], the integration kernel of eq. ()52|) can be expanded in 
a Fourier series on the group 

exp (^[W + TtV^] + ^TtV Trl-t j = ^ c,x^^^\V) , (53) 

where the sum runs over the set of all irreducible representations of SU(3) 
and x^'^\^) is the character of the representation u. In order that eq. 
holds, it is necessary and sufficient that all coefficients > Q are positive. 
For /?/ > and /3a = it is proven that this is indeed the case: 

e^v(^[W + TiV^]] = f; a„^(Tr\/)"(Tr\/t)- ^ (54) 



n,m=0 



where a„m > 0. Here, (Try)"(Try^)™ is the trace of the tensor product rep- 
resentation of SU(3) composed of n quark and m antiquark representations. 
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Reducing out the tensor product, one gets 



(Tr\/)"(Tr0)™ = J2cu{n,m)x^''\V) , (55) 

V 

where Ci^{n,m) > 0. Since all irreducible representations can be obtained by 
reducing out tensor products of quark representations, > for all z/. It is 
trivial to extend this prove to /9/ > and /3a > 0. 

However for (3a < 0, as we consider here, the situation becomes more 
complicated. In the expansion 

exp(^^[Try + Wt] + ^TrV TrVt^ = E a„„( TrV)'^( W^)™ (55) 

it is no longer guaranteed that a„m > for all choices of n, m. This can be 
most easily seen for 



2N 



n+l 



1 

H — tU 
nl 



2N 



n-l 



(3a_ 

ATS 



(57) 



It follows that ttn,! is positive for (/?j/2)^/n + /9a > 0. I.e. as n increases, the 
lower limit on j3a goes to zero. However, it remains quite unclear how the 
(in,m will add up in the coefficients Ci, and whether > 0. We were not able 
to clarify this question rigorously. 

To get some idea, we evaluated the coefficients for the pairs of Pf,Pa 
studied in this paper, for representations u up to the dimension 15. 

To this end we evaluated the integrals 

Cu = J dVxtiV) exp (^^[Try + Tr0] + ^TiV TrV^^ (58) 

numerically. It turned out that > for all coefficients that we computed, 
except for the coefficient of one 15 dimensional representation at f3f = 7.8477 
and Pa = —4. 

We applied a second numerical approach to check the positivity of the 
integration kernel. Assume that we perform the integration eq. with 
a Monte Carlo method. I.e. we evaluate f{U) for m SU(3) matrices that 
have been selected randomly. Choosing the same SU(3) matrices for both 
integrations, the integration kernel becomes a real symmetric mx m matrix. 
In this study, we have used m < 10000. 

To check the positivity of this matrix, we evaluated its smallest eigenvalue. 
It turned out that for = —2 and Pa = —4 and all values of Pf that we 
have studied here, negative eigenvalues are found. However, in particular for 
Pa = —2, the absolute value of the smallest eigenvalue is by several orders of 
magnitude smaller than the largest eigenvalue. Notice that no quantitative 
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sign of a violation of positivity has been observed in the decay of the Wilson 
loop correlation functions, while for the improved actions studied in ref. jl] 
such violations were clearly visible. It remains an open question, whether 
there is a finite range of negative jSa, where the transfermatrix is strictly 
positive. 
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